Kinetic step bunching instability during surface growth 
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We study the step bunching kinetic instability in a growing crystal surface characterized by 
anisotropic diffusion. The instability is due to the interplay between the elastic interactions and 
the alternation of step parameters. This instability is predicted to occur on a vicinal semiconductor 
surface Si(001) or Ge(001) during epitaxial growth. The maximal growth rate of the step bunching 
increases like F 4 , where F is the deposition flux. Our results are complemented with numerical 
simulations which reveals a coarsening behavior on the long time for the nonlinear step dynamics. 
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The field of surface growth on semiconductors is very 
active due to its importance for both technological ap- 
plications and fundamental science [H 0- 13 Under non- 
equilibrium growth a variety of experiments reveal rich 
crystal morphologies resultingfrom the nonlinear evolu- 
tion of step flow instabilities jj, 0| . The kinetic and elas- 
tic effects drive the system towards self-organized states, 
characterized by the appearance of ordered structures on 
the vicinal surfaces. This self-organization can be ex- 
ploited in semiconductor nanotechnology, for the devel- 
opment of devices having interesting quantum properties 
[a, ■ A basic mechanism for pattern formation in vici- 
nal semiconductor surfaces is the step bunching instabil- 
ity 0- 0- E3 whose origin is commonly attributed to the 
presence of impurities, to the inverse Ehrlich-Schwoebel 
effect, or to electro- migration (see [Tllllij and references 
therein) . This Letter is motivated by the molecular beam 
epitaxy (MBE) experiments on Si(001) in which it was 
observed a new type of kinetic instability leading to the 
formation of step bunches [151 H3. flij | . Microscopic ki- 
netic Monte-Carlo simulation showed that in the case of 
Si(001), step bunching is due to the coupling between 
diffusion anisotropy and differences in step kinetic pa- 
rameters . Here we provide a macroscopic instability 
mechanism for the step bunching instability that does not 
require any inverse Ehrlich-Schwoebel effect. We show 
that the interplay between the elastic step interactions 
and the alternation of kinetic parameters, characteristic 
of the Si(001) vicinal surface, induces a finite wavelength 
instability with maximal growth rate increasing as F 4 (F 
is the deposition flux) . Our results are complemented by 
numerical simulations which reveal a coarsening behavior 
on the long time for the non-linear step dynamics. 

The Si(001) vicinal surface consists of a periodic se- 
quence of terraces where rows of 2 x 1 dimerised adatoms 
(terrace of type a) alternate with 1x2 dimerised adatoms 
(terrace of type 6), see Fig.Q] On the reconstructed sur- 
face adatoms diffuse preferentially along dimcr rows, giv- 
ing rise to an anisotropic diffusion. Therefore, the steps 
separating the terraces are of two kinds. The S a step is 
rather straight while the Sb step is very corrugated 




FIG. 1: Sketch of the Si(001) vicinal surface showing the al- 
ternation of terraces and steps S a and St- D a and D(, are the 
surface diffusion coefficients, and v the upper and lower ter- 
race kinetic attachment coefficient. The horizontal axis shows 
terrace lengths 1% and step positions x™, a;™. 



The significant difference in the attachment kinetics be- 
tween smooth and rough steps allows us to choose the 
kinetic coefficient v a — v of the S a step to be finite and 
the one of the Sb step to be infinite [yb — > oo) |l8j . 

Let us denote by x™(t) and x^(t) the positions at time 
t of steps S a and Sb respectively, the terrace sizes are 
'?(*) = - and l%{t) = a£(i) - x»(t), and 

v2(t) — i"(t), v£(t) = x™(t) are the step velocities, where 
n = 1,2,... numbers the terrace (Fig. 0. During the 
growth, the adatom concentration on each terrace obeys 
the following diffusion equations [l9|. 



D "^2 C a(*^) = -Fi and D b -^C^(x,t) = -F, (1) 



where C™(x,t) and C^(x,t) are the adatom concentra- 
tion on each terrace, D a and Db the diffusion coefficients, 
and F the deposition flux. We assume that the des- 
orption of adatoms is negligible, since step bunching is 
observed at temperatures less then 1200 K; we also ne- 
glect transparency of steps. Equations are subject to 
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FIG. 2: Terrace size h as a function of the flux / (nondimen- 
sional) for the fixed point solution with a a = 1 and ctb = 0.1. 
For / > f c = 4(2 + a a + a b )(3 + 2a a - 2 ^(1 + a a )(2 + a a ) ) 
no stationary solution exists. For f < fc = 4a&/(l + a a +a!&) 
the upper branch disappears. 



boundary conditions at each step: 

(~m( n .n— 1\ _ /-ro 



.arm 

9cr 



<9x 



«) = ^[C 6 "«)-C, 



(2) 
(3) 

(4) 
(5) 



where v is the step kinetic coefficient of step a; C™ g a and 
C" q b are the adatom equilibrium concentration at each 
step. The velocity of each step is given by the conditions, 
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dx dx 
with f2 the unit atomic surface. The adatom equilibrium 
concentrations are determined by the elastic interactions 
between steps mediated by the terrace reconstruction [Til 

MM, 



eq,b 
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where E = ^ICoA/hsT (Co is the uniform equilibrium 
concentration, ks the Boltzmann constant, T the tem- 
perature, and A is an energy per unit length related to the 
stress and elastic constants of the medium 213). The sys- 
tem behavior is controlled by three independent nondi- 
mensional parameters: / = Fl^/vE, a a = vl a /D a and 
otb = vla/Db- We set the unit length to be the initial 
size of terrace a, l a = 1. With this convention the ve- 
locity is measured in units of vECl/l a . Under normal ex- 
perimental conditions (£IF < Is -1 , A w 5 1Q — 12 Jm — 1 , 
nC w 1CT 2 ) / < 1, a a > 1 and a b < 1. 
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FIG. 3: The growth rate s versus k, (a) real and (b) imaginary 
parts, for the 1^ branch with a a = 1, at — 0.1, and / = 1.2 < 

U 



We first look for a uniform train of steps with Z™ = 1 
and l b — lb for all n, traveling at a constant velocity. The 
relative velocity is 

_ 4-4(4-/) , 4-J 6 (4-Z 6 /) 



l b (l + a a ) l b (l + l b a b ) 
It vanishes for terraces of size lb giving by 
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where A = (2 - //4 + a a - a b ) 2 - /(2 + a«)(l + a a + 
ab — 4ctb/f). Depending on the values of the parameters 
(a a ,ctb,f) different situations may arise. The two solu- 
tions l^ represent terraces of nearly equal size {l b ), or 
double steps (l b ) for which terraces of type a almost dis- 
appear (Fig. |2l . This two branches exist (A > 0, lb > 0) 
for fluxes lower than f c . The upper branch characteriz- 
ing the double step state, exists in the interval (/^,/ c ), 
as shown in Fig. |2 This interval vanishes for the values 
of ab satisfying a\ > (1 + a a )(2 + a a ). 

We investigate now the stability of a train of steps 
traveling at the constant velocity V = v a (lb) = Vb(h), 
where lb is one of the solutions (|11|) . The perturbed step 
positions are given by = Vt + n(l + l b ) + X™ and 
x b = Vt + n (l + h) + h + x b- We expand the step 
velocities up to first order in the perturbation X£ b (t) = 
X a (,exp(st -+- ink), where s is the (complex) growth rate 
and k the wavenumber. After some algebra, we obtain 
the following expression for the growth rate s = a + iw, 
in the limit of small / and k, for the lower branch l b , 
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FIG. 4: Maximum growth rate 5Rs = a ma x versus /, for the 
branch with ot a = 1, and ot b = 0.1. Inset in log-log scale 
showing the / 4 law for small /. 

where 

1 (q + a b + 3) 
16(a Q + a b + 2) 2 



1536(a a + a h + 2) 3 
[768(1 + a Q + a h ) + 192(a a + a b ) 2 - 
24(8 + 6a a + Ua b + a 2 a + 6a a a b + 5a 2 ) f - 
(72 + 86(a a + at) + I9a 2 a + 80a a a h - Hag)/ 2 ] 
_ 1 8a a + 8a b + 16 - 9(a a - a b )/ 
' _ 384 2 + a a + a fc 

Figure [21 shows the real and imaginary parts of s for 
typical values of the parameters. The imaginary part 
lu = 3s ~ fk is associated with the propagation of com- 
pression waves. The upper branch l£ displays a similar 
behavior, showing that both uniform solutions are unsta- 
ble. The growth rate a = 5is of this instability increases 
like f 2 k 2 for small k and is stabilized at larger k by elas- 
tic effects (independent of /). The most unstable mode 
kmax increases linearly with the flux k max ~ /, and its 
corresponding growth rate scales like a max ~ / 4 , as can 
be seen from Eq. (|12fl and illustrated in Fig. 21 Therefore, 
for small /, the instability is almost vanishing. We have 
checked that this result is not qualitatively changed if a 
finite kinetic coefficient v b for the b step is included into 
the model. In the particular case v a = v b , the growth 
of the step flow instability vanishes, independently of the 
D a /D b ratio. 

This step bunching instability is a consequence of ki- 
netic and elastic coupling under crystal growth condi- 
tions. The asymmetry of kinetic coefficients, related to 
the alternation of smooth and rough steps, introduces a 
supplementary elastic coupling between non-neighboring 
terraces. Their combined action breaks the symmetry 
between the ascending and descending adatoms currents. 




FIG. 5: Initial spacetime evolution of a train of iV = 40 steps 
x„(t) (mean velocity subtracted). Numerical resolution of 
Eqs. 1131141 with a a = 1, a b = 0.1, and / = 1.2. 
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FIG. 6: Snapshots of the surface during step bunching for 
a a = 1 and a b = 0.1, / = 1.2 and N = 128. (a) The mod- 
ulated step pattern at t — 400; (b) the coarsening regime at 
t = 1000. 



In contrast to the usual inverse Ehrlich-Schwoebel effect, 
the unbalanced currents concern the interaction between 
two double terraces. 

We can obtain the non-linear evolution of the step po- 
sitions by numerical integration of the set of differential 
equations I|6I7|) . These equations can be written explic- 
itly using Q and its associated boundaries conditions 
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iPEJ, they read: 



1 + «b^ 1 
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where M 2 = l/l2 +1 +l/J?-2/l?, and //" = l/ij^+l/i?- 
2/7™. The typical evolution of the step flow on the vicinal 
surface is shown in Figs. [5] and |S] We choose as initial 
conditions different sets of slightly perturbed equal size 
terraces. The amplitude of the perturbations increases 
during the evolution while compression waves propagate 
as predicted by the linear analysis. The further evolu- 
tion of the system shows the formation of bunches of 
steps separated by wider terraces. A coarsening regime 
sets in, characterized by the coalescence of these step 
bunches into larger ones. Ultimately the long time evo- 
lution leads to a surface composed of a single macro-step. 
We have checked that these results are not affected by the 
choice of initial conditions or by the size of the system. 
This behavior is illustrated in Fig. by numerical sim- 
ulations depicting the step train at two different times 
for f < f c - We have also checked that this scenario re- 
mains valid when / > f c . FigureUJa) shows a modulated 
step pattern characteristic of the incipient formation of 
step bunches. Figure EJb) shows the further evolution of 
these bunches during the coarsening regime. A variety 
of models of step bunching shows that the mean number 
of step s in a single bunch obeys a power law (n(t)) ~ 
pj-lllj. In the present case this behavior is expected to be 
more complex due to the presence of oscillations in the 
step size. This point deserves a specific investigation us- 
ing a continuous description of the step dynamics, more 
appropriated for the study of the long time evolution. 

In conclusion, we have proposed a kinetic mechanism of 
step bunching instability under non-equilibrium growth 
of a vicinal surface of Si(001)-like crystals. The essen- 
tial ingredients are the alternation of rough and smooth 
steps, having different diffusion and attachment coeffi- 
cients, and the coupling of non-neighboring terraces by 
elastic interactions. We have found that the rate of 
growth of the instability scales like / 4 . It would be inter- 
esting to experimentally test this behavior, using for ex- 
ample measurements of the surface roughness evolution. 
Some indirect evidence supporting the present instability 
mechanism are provided by experiments on Si/Ge growth 
[2^| . The addition of a small quantity of Ge atoms pro- 
duces a substantial modification of the steps properties, 
and in particular shades off the distinction between rough 
and smooth steps. As a consequence, a stable step flow 
regime sets in even though the diffusion coefficient of 
successive terraces are still different. A generalization 



of the present model to two dimensions is currently un- 
der progress. We expect a rich phenomenology due to 
the coupling between a meandering instability induced 
by the alternation of the step parameters [13 , the step 
bunching instability and the anisotropy of surface diffu- 
sion H3. 
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